#ifndef KERNEL_PROBIT_HPP
#define KERNEL_PROBIT_HPP
#include "../Random.hpp"
#include "../simulate/rtnorm.hpp"
#include "../simulate/rtnomNC.hpp"
#include <armadillo>

using namespace arma; 

template<class Density>//Proposal avec r et d
class Kernel_Probit{
	public:
		Kernel_Probit(Density *d)
		{
			_d=d;
			//_p=malloc(sizeof(Proposal));
			boost::random::uniform_01<> U;
			runif=new RandomG::Random<boost::random::uniform_01<> >(U);
			rt=new rtnormNC();
			_Sig=1;
			_X=d->Get_X();
			_Y=d->Get_Y();
			_p=_X.n_cols;
			_n=_X.n_rows;
			_v=d->Get_v();
		}
		~Kernel_Probit(){
			delete runif;
			delete rt;
		}
/*		Kernel_MH<Density,Proposal>& Kernel_MH<Density,Proposal>::operator=(const Kernel_MH<Density,Proposal>& X){
			cout << "Hello";	
		}
*/
		void Move(mat *X,int i)//Move matrix convert to multidimentional theta 
		{		//, moves each column dimentionis in line
			double xi=_d->Get_bn();	
			int n=(*X).n_rows;
			mat u(_n,_n);
			u.fill(pow(xi,2));
			_V=diagmat(u);	
			for(int i=0;i<n;i++)
			{
				mat z=rz(((*X)(i,span::all)).t(),xi);
				mat b=Pi_b(z);
				(*X)(i,span::all)=b.t();
				//cout << i << "\n";
			}
			//////////
			//delete[] y;
		}

		mat Pi_b(mat z)
		{
			mat temp=_X.t()*inv(_V)*_X;
			mat V=inv(inv(_v)+temp);
			mat m=V*_X.t()*inv(_V)*z;	
			/*cout << V;
			cout << X2;
			cout << nb;*/
			G=new Distribution::Gaussian(_p,m,V);
			mat res1=(*G).r(1);
			//cout << z.t();
			//cout << beta.t();	
			delete G;	
			return res1;			
		}	
		mat rz(mat b,double xi)
		{
			mat z(_n,1);
			for(int i=0;i<_n;i++)
			{
				double foo=dot(_X.row(i),b);
				if(_Y[i]==0)
				{
					z(i,0)=(*rt)(foo,-pINF,0,xi);
					//cout << z(i,0) << " y " << _Y[i] << "\n"; 
				}else{
					z(i,0)=(*rt)(foo,0,pINF,xi);
				}
			}
			return z;
		}

		void Set_bool(int b){_bool=b;}
		void Set_n(int b){
		}
		
	
	private:
		int _sum;
		Density *_d;
		int _bool;
		rtnormNC *rt;
		mat _X;
		double *_Y;
		mat _v;
		int _p;
		int _n;
		Distribution::Gaussian *G;
		RandomG::Random<boost::random::uniform_01<> > *runif;
		double _Sig;
		int _i;
		mat _V;
};
#endif

